function idat = fsb_enhance_contrast(idat,sandbox)

% FSB : improve contrast between gray and white matter
%
% EXAMPLE:
% idat = fsb_enhance_contrast2(idat,sandbox)
%
% INPUT:
% idat: 4 D image data
% sandbox: sandbox experiment struct
%
% OUTPUT:
% idat: contrast enhanced 4D image data
%
% CALLED BY:
% FSB.m
%
% NOTES:
% Might be useful if you have to coregister different modalities
%
% Copyright 2010 MPI for Biological Cybernetics
% Author: Steffen Stoewer
% License:GNU GPL, no express or implied warranties
% 
% $Revision 1.0
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% Check for size of dataset
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

d = size(idat);
if size(d)<4;
    disp('No 4D dataset');
    sandbox=[];
end

idat = single(idat);

%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% Check if sandbox struct exists
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

if isempty(sandbox)

    %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    % Do enhancement for single volume
    %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    disp('Enhancing contrast ...')
    h = waitbar(0,'Enhancing contrast...');
    for x = 1:d(1);
        waitbar(x/d(1));
        for y = 1:d(2);
            te = single(squeeze(idat(x,y,:,:)));
            meante = repmat(mean(te,2),[1,size(te,2)]);

            det = int16(detrend(te)+meante);
            idat(x,y,:,:)=det;
        end
    end
    close(h);
    disp('...done');
    hold on;

else

    %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    % Go over runs
    %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    for scannum = 1:max(sandbox.intrial(:,4));

        detstring = ['Enhancing contrast in run ' num2str(scannum) ' of ' num2str(max(sandbox.intrial(:,4)))];
        disp(detstring)

        h = waitbar(0,detstring);
        scanindex = find(sandbox.intrial(:,4) == scannum);
        idat_seg = idat(:,:,:,scanindex);

        %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
        % Go over volumes
        %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
        for x = 1:d(1);
            waitbar(x/d(1));
            for y = 1:d(2);

                te = single(squeeze(idat_seg(x,y,:,:)));
                meante = repmat(mean(te,2),[1,size(te,2)]);
                
                %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
                % Enhance contrast
                %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
                det = int16(detrend(te)+meante);
                idat_seg(x,y,:,:)=det;
            end
        end

        close(h);
        disp(['...done with run ' num2str(scannum)]);
        idat(:,:,:,scanindex) = idat_seg;

    end

end

end
